Binary black hole merger in the extreme- mass-ratio limit: a multipolar analysis. 
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Building up on previous work, we present a new calculation of the gravitational wave emission 
generated during the transition from quasi-circular inspiral to plunge, merger and ringdown by a 
binary system of nonspinning black holes, of masses mi and m2, in the extreme mass ratio limit, 
mim2 <^ {mi + m^)^ . The relative dynamics of the system is computed without making any adia- 
batic approximation by using an effective one body (EOB) description, namely by representing the 
binary by an effective particle of mass jj, — m\m2/(mi + 7712) moving in a (quasi- )Schwarzschild 
background of mass M = mi -I- m2 and submitted to an 0{v) 5PN-resummed analytical radia- 
tion reaction force, with v = [ijM. The gravitational wave emission is calculated via a multipolar 
Regge- Wheeler- Zerilli type perturbative approach (valid in the limit <C 1). We consider three 
mass ratios, v = {10~^, 10"'', 10^*}, and we compute the multipolar waveform up to ^ = 8. We esti- 
mate energy and angular momentum losses during the quasiuniversal and quasigeodesic part of the 
plunge phase and we analyze the structure of the ringdown. We calculate the gravitational recoil, 
or "kick", imparted to the merger remnant by the gravitational wave emission and we emphasize 
the importance of higher multipoles to get a final value of the recoil v/{cv^) = 0.0446. We finally 
show that there is an excellent fractional agreement (~ 10""^) (even during the plunge) between 
the 5PN EOB analytically-resummed radiation reaction fiux and the numerically computed gravi- 
tational wave angular momentum flux. This is a further confirmation of the aptitude of the EOB 
formalism to accurately model extreme-mass-ratio inspirals, as needed for the future space-based 
LISA gravitational wave detector. 

PACS numbers: 04.30. Db, 04.25. Nx, 95.30.Sf, 



I. INTRODUCTION 

After the breakthroughs of 2005 [1-3] , state-of-the-art 
numerical relativity (NR) codes can nowadays routinely 
evolve (spinning) coalescing binary black hole systems 
with comparable masses and extract the gravitational 
wave (GW) signal with high accuracy [4-11]. However, 
despite these considerable improvements, the numerical 
computation of coalescing black hole binaries where the 
mass ratio is considerably different from 1:1 is still chal- 
lenging. To date, the mass ratio 10:1 (without spin) re- 
mains the highest that was possible to numerically evolve 
through the transition from inspiral to plunge and merger 
with reasonable accuracy [12, 13]. In recent years, work 
at the interface between analytical and numerical rel- 
ativity, notably using the effective-one-body (EOB) re- 
summed analytical formalism [14-19], has demonstrated 
the possibility of using NR results to develop accurate an- 
alytical models of dynamics and waveforms of coalescing 
black-hole binaries [20-28]. 

By contrast, when the mass ratio is large, approxima- 
tion methods based on black-hole perturbation theory 
are expected to yield accurate results, therefore enlarg- 
ing our black-hole binaries knowledge by a complemen- 
tary perspective. In addition, when the larger black-hole 
mass is in the range 10^Mq-10''Mq, the GWs emitted 
by the radiative inspiral of the small object fall within 
the sensitivity band of the proposed space-based detec- 
tor LISA [29, 30], so that an accurate modelization of 
these extreme-mass-ratio-inspirals (EMRI) is one of the 
goals of current GW research. 



The first calculation of the complete gravitational 
waveform emitted during the transition from inspiral to 
merger in the extreme-mass-ratio limit was performed 
in Refs. [31, 32], thanks to the combination of 2.5PN 
Fade resummed radiation reaction force [33] with Regge- 
Wheeler-Zerilli perturbation theory [34-37]. This test- 
mass laboratory was then used to understand, element 
by element, the physics that enters in the dynamics and 
waveforms during the transition from inspiral to plunge 
(followed by merger and ringdown) , providing important 
inputs for EOB-based analytical models. In particular, 
it helped to: (i) discriminate between two expressions 
of the resummed radiation reaction force; (ii) quantify 
the accuracy of the resummed multipolar wavefom; (iii) 
quantify the effect of non-quasi-circular corrections (both 
to waveform and radiation reaction); (iv) qualitatively 
understand the process of generation of quasi-normal 
modes (QNMs); and (v) improve the matching procedure 
of the "insplunge" waveform to a "ringdown" waveform 
with (several) QNMs. In the same spirit, the multipo- 
lar expansion of the gravitational wave luminosity of a 
test-particle in circular orbits on a Schwarzschild back- 
ground [38-40] was helpful to devise an improved resum- 
mation procedure [21, 41] of the PN (Taylor-expanded) 
multipolar waveform [42, 43]. Such resummation proce- 
dure is one of the cardinal elements of what we think is 
presently the best EOB analytical model [26, 27]. Sim- 
ilarly, Ref. [44] compared "calibrated" EOB-resummcd 
waveforms [24] with Teukolsky-based perturbative wave- 
forms, and confirmed that the EOB framework is well 
suited to model EMRIs for LISA. 
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In addition, recent numerical achievement in the calcu- 
lation of the conservative gravitational self-force (GSF) 
of circular orbits in a Schwarzschild background [45-47] 
prompted the interplay between post-Newtonian (PN) 
and GSF efforts [48, 49], and BOB and GSF efforts [50]. 
In particular, the information coming from GSF data 
helped to break the degeneracy (among some EOB pa- 
rameters) which was left after using comparable-mass 
NR data to constrain the EOB formalism [50]. (See also 
Ref. [28] for a different way to incorporate GSF results 
in FOB). 

In this paper we present a revisited computation of 
the GWs emission from the transition from inspiral to 
plunge in the test-mass limit. We improve the previous 
calculation of Nagar et al. [31] in two aspects: one nu- 
merical and the other analytical. The first is that we use 
a more accurate (4th-order) numerical algorithm to solve 
the Regge- Wheeler- Zerilli equations numerically; this al- 
lows us to capture the higher order multipolar informa- 
tion (up to £ = 8) more accurately than in [31]. The 
second aspect is that wc have replaced the 2.5PN Fade 
resummed radiation reaction force of [31] with the 5PN 
resummed one that relies on the results of Ref. [41]. 

The aim of this paper is then two-fold. On the one 
hand, our new test-mass perturbative allows us to de- 
scribe in full, and with high accuracy, the properties of 
the gravitational radiation emitted during the transition 
inspiral-plunge-merger and ringdown, without making the 
adiabatic approximation which is the hallmark of most 
existing approaches to the GW emission by EMRI sys- 
tems [51, 52]. We compute the multipolar waveform up 
to € = 8, we discuss the relative weight of each multipole 
during the nonadiabatic plunge phase, and describe the 
structure of the ringdown. In addition, from the multi- 
polar waveform we compute also the total recoil, or kick, 
imparted to the system by the wave emission, thereby 
complementing NR results [12]. On the other hand, we 
can use our upgraded test-mass laboratory to provide 
inputs for the EOB formalism, notably for completing 
the EOB multipolar waveform during the late-inspiral, 
plunge and merger. As a first step in this direction, we 
show that the analytically resummed radiation reaction 
introduced in [26, 41] gives an excellent fractional agree- 
ment (~ 10~^) with the angular momentum flux com- 
puted a la Regge- Wheeler- Zerilli even during the plunge 
phase. 

This paper is organized as follows. In Sec. II we give a 
summary of the formalism employed. In Sec. Ill we de- 
scribe the multipolar structure of the waveforms; details 
on the energy and angular momentum emitted during 
the plunge-merge-ringdown transition are presented as 
well as an analysis of the ringdown phase. Section IV 
is devoted to present the computation of the final kick, 
emphasizing the importance of high multipoles. The fol- 
lowing Sec. V is devoted to some consistency checks: on 
the one hand, we discuss the aforementioned agreement 
between the mechanical angular momentum loss and GW 
energy flux during the plunge; on the other hand, we in- 



vestigate the influence of EOB "self-force" terms (either 
in the conservative and non conservative part of the dy- 
namics) on the waveforms. We present a summary of our 
findings in Sec. VI. In Appendix A we supply some tech- 
nical details related to our numerical framework, while 
in Appendix B wc list some useful numbers. We use ge- 
ometric units with c = G ~ 1 . 



II. ANALYTIC FRAMEWORK 
A. Relative dynamics 

The relative dynamics of the system is modeled spec- 
ifying the EOB dynamics to the small-mass limit. The 
formalism that we use here is the specialization to the 
test-mass limit of the improved EOB formalism intro- 
duced in Ref. [26] that crucially relies on the "improved 
resummation" procedure of the multipolar waveform of 
Ref. [41] . Let us recall that the EOB approach to the gen- 
eral relativistic two-body dynamics is a nonperturhatively 
resummed analytic technique which has been developed 
in Refs. [14-18]. This technique uses, as basic input, the 
results of PN theory, such as: (i) PN-cxpanded equations 
of motion for two pointlike bodies, (ii) PN-expanded ra- 
diative multipole moments, and (iii) PN-expanded en- 
ergy and angular momentum fluxes at infinity. For the 
moment, the most accurate such results are the 3PN con- 
servative dynamics [53, 54], the 3.5PN energy fiux [55- 
57] for the v ^ Q case, and 5.5PN [58] accuracy for the 
V = Q case. Then the EOB approach "packages" this PN- 
expanded information in special resummed forms which 
extend the validity of the PN results beyond the expected 
weak-field-slow-velocity regime into (part of) the strong- 
field-fast-motion regime. In the EOB approach the rela- 
tive dynamics of a binary system of masses mi and m2 is 
described by a Hamiltonian H-eob{M^ /i) and a radiation 
reaction force J^eob(-W, m), where M = mi+m2 and fi = 
mim2/M. In the general comparable-mass case Heob 

has the structure HeobI-^-'^, m) = A/y^ 1 + 2v{H^ — 1) 
where v = fi/M = mim2 / {mi + m2)^ is the symmetric 
mass ratio. In the test mass limit that we are considering, 
v ^ 1, we can expand Heob in powers of v. After sub- 
tracting inessential constants we get a Hamiltonian per 
unit (^) mass H = lim^^o(-ff — const. )/^ = lim^^Q H^. 
As in Refs. [31, 32], we replace the Schwarzschild radial 
coordinated = r + 2Mlog[r/(2M) — 1] and, corre- 
spondingly, the radial momentum by the conjugate 
momentum Pr^ of i?*, so that the specific Hamiltonian 
has the form 



h = ]]a[i + ^^^+pI . (1) 

Here wc have introduced dimensionless variables ?" = 
R/M,f^ = i?,/Af , pr, = PrJh, p^ = Pv/ifJ-M) and A = 
1 — 2/?". Hamilton's canonical equations for (f, r, , p^, , P<p) 
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in the equatorial plane {9 = tt/2) yield 
: Pf, 

A A 

r = —Pr, = Vr , 

H 

^Pv -n 



Pr, 



r-2 



(2) 

(3) 

(4) 

(5) 
(6) 



Note that the quantity is dimensionless and represents 
the orbital frequency in units of 1/M. In these equations 
the extra term [of order Oiy)] represents the non 
conservative part of the dynamics, namely the radiation 
reaction force. Following [31, 32, 73], we use the following 
expression: 
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(7) 



where Vtp — ffi is the azimuthal velocity and /din = 
/Fncwi denotes the (Newton normalized) energy 
flux up to multipolar order ^max(in the v = Q limit) re- 
summed according to the "improved resummation" tech- 
nique of Ref. [41]. This resummation procedure is based 
on a particular multiplicative decomposition of the multi- 
polar gravitational waveform. The energy flux is written 
as 



1 

8^ 



1=2 m=l 



where hpm is the factorized waveform of [41] 



I'lm 



Tim 6 



Per, 



(8) 



(9) 



where h^g^''' represents the Newtonian contribution given 
by Eq. (4) of [41], e = (or 1) for ^ + m even (odd), 
source", Eqs. (15-16) of [41]; Tf„ 



Sgjj is the effective 



is the "tail factor" that resums an infinite number of 
"leading logarithms" due to tail effects, Eq. (19) of [41]; 
6i,n is a residual phase correction, Eqs. (20-28) of [41]; 
and pim is the residual modulus correction, Eqs. (Cl- 
C35) in [41]. In our setup we truncate the sum on £ at 
•^max = 8. We refer the reader to Fig. 1 (b) of [41] to figure 
out the capability of the new resummation procedure to 
reproduce the actual flux (computed numerically) for the 
sequence of circular orbits in Scliwarzshild^. 



^ We mention that, although Ref. [41] also proposes to further 
(Pade) resum the residual amplitude corrections (and, in 

particular, the dominant one, P22) to improve the agreement 
with the "exact" data, we prefer not to include any of these so- 



B. Gravitational wave generation 

The computation of the gravitational waves gener- 
ated by the relative dynamics follows th same line of 
Refs. [31, 32], and relies on the numerical solution, in 
the time domain, of the Regge- Wheeler- Zerilli equations 
for metric perturbations of the Schwarzschild black hole 
with a point-particle source. Once the dynamics from 
Eqs. (2)-(6) is computed, one needs to solve numerically 
(for each multipole {£, m) of even (e) or odd (o) type) a 
couple of decoupled partial differential equations 



(c/o) 



(0/0) 



(10) 



with source terms S^^J^°^ linked to the dynamics of the 
binary. Following [31], the sources are written in the 
functional form 



+ Ft':\r,t)dr.5{r*-R*{t)) 



(11) 



with r-dependent [rather than i?(i)-dependent] coeffi- 
cients G{r) and F{r). The explicit expression of the 
sources is given in Eqs. (20-21) of [31], to which we ad- 
dress the reader for further technical details. We men- 
tion, however, that in our approach the distributional 
i5-function is approximated by a narrow Gaussian of fi- 
nite width a ^ M . In Ref. [31] it was already pointed 
out that, if a is sufficiently small and the resolution is 
sufficiently high (so that the Gaussian can be cleanly 
resolved) this approximation is competitive with other 
approaches that employ a mathematically more rigorous 
treatment of the (5-function [67-69] (see in this respect 
Table 1 and Fig. 2 of Ref . [31]). That analysis motivates 
us to use the same representation of the J-function also 
in this paper, but together with an improved numerical 
algorithm to solve the wave equations. In fact, the solu- 
tion of Eqs. (10) is now provided via the method of lines 
by means of a 4th-order Runge-Kutta algorithm with 
4th-order finite differences used to approximate the space 
derivatives. This yields better accuracy in the waveforms 
(using resolutions comparable to those of Ref. [31]), and 
allows to better resolve the higher multipoles. More de- 
tails about the numerical implementation, convergence 
properties, accuracy, and comparison with published re- 
sults are given in Appendix A. 

From the numerically calculated master functions 



(c/o) 



one can then obtain, when considering the limit 



phistications here. This is motivated by the fact that, along the 
sequence of circular orbits, all the different choices are practically 
equivalent up to (and sometimes below) the adiabatic last stable 
orbit (LSO) at r = 6M (see in this respect their Fig. 5). In prac- 
tice our P22 actually corresponds to the Taylor-expanded version 
(at 5PN order) of the remnant amplitude correction, denoted 
n[p22] in [41]. 
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r — (X), the /i+ and hx gravitational- wave polarization 
amplitude 



h+ — i/i> 



£>2,rn 



(^-2)! 



Ini] -2J^ 



(12) 



where -zF^™ = -2F^™(6',(p) are the s = 2 spin- 
weighted spherical harmonics [59] . From this expression, 
all the interesting second-order quantities follow. The 
emitted power, 



(^ + 2)! 



IGtt ^ (^-2)! 



the angular momentum flux 



^(o)^(o)* , ^(o)^(o)* 



and the linear momentum flux [70-72] 



(13) 

-t-c.c. 
(14) 



£>2,m 



^(0)4,(0)* 



(15) 



with 



a^™ = 2(^ - 1)(£ + 2)v/(^ - m)(^ + m + 1) (16) 
(^-t-3)! /(^ + m-|- l)(£-hm + 2) 



3£m 



(£+l)(£-2)! 



(2£+ l)(2£-h3) 



(17) 



III. RELATIVE DYNAMICS AND WAVEFORMS 

Let us now consider the dynamics and waveforms ob- 
tained within our new setup. Evidently, at the qualitative 
level our results are analogous to those of Rcfs. [31, 32]. 
By contrast, at the quantitative level, dynamics and 
waveforms are slightly different due to the new, more 
accurate, radiation reaction force. The particle is ini- 
tially at r = 7M. The dynamics is initiated with the so 
called post-circular initial data for (ptpjPr) introduced in 
Ref. [15] and specialized to the fi —> limit (see Eqs. (9)- 
(13) of Rcf. [31]). Because of the smallness of the value of 
yU we are using, this approximation is sufficient to guar- 
antee that the initial eccentricity is negligible. To have 
a better modelization of the extreme-mass-ratio limit 
regime we considered three values of the mass ratio ly, 
namely v = {lO'^, 10"^, 10"''}. The values of v are 
chosen so that the particle passes through a long (when 
v < 10~^) quasicircular adiabatic inspiral before enter- 
ing the nonadiabatic plunge phase. Fig 1 displays the 
relative trajectory for ly = 10"'^. The system executes 




FIG. 1. Transition from quasicircular inspiral orbit to plunge. 
Initial position is ro = 7M and v — 10~'^. 



about 40 orbits before crossing the LSO at r = 6M while 
plunging into the black hole. 

The main multipolar contribution to the gravitational 
signal is clearly the £ = m = 2. The real part of the corre- 
sponding waveform is displayed in Fig. 2. It is extracted 



at r' 



obs 



lOOOA'/ and it is shown versus observer's re- 



tarded time u = {t°^^ — r°^^)M . Note how the amplitude 
of the long wavetrain emitted during the adiabatic quasi- 
circular inspiral grows very slowly for about 4000il/, until 
the transition from inspiral to plunge around the crossing 
of the adiabatic LSO frequency. In the following we want 
however to focus on the higher order multipolar contribu- 
tions to the waveform, as they are particularly relevant 
in our test-mass setup. The computation of these multi- 
poles and their inclusion in the analysis is one of the new 
results of this paper ^. 

Figure 3 summarizes the information a bout the mul- 
tipolar waveforms up to ^ = 8. The left panels show 
the moduli (normalized by the mass ratio v), while the 
right panels show the corresponding instantaneous gravi- 
tational wave frequencies MuJem- We show, for each value 
of £, the dominant (even-parity) ones, i.e. those with 
m — i, together with some subdominant (odd-parity) 
ones. The comparison between the moduli highlights how 
the amplitude of higher modes, that is almost negligible 
during the adiabatic inspiral, can be magnified of about 
factor two (see the £ = 2, m = 1 case) or three (see the 
m = £ = 6 case) during the nonadiabatic plunge phase. 
This fact is expected to have some relevance in those 
computations that are dominated by the nonadiabatic 



^ We note that calculations up to £ = 4 were already performed 
in Ref. [31, 32], but no higher-order multipolar waveforms were 
either shown or discussed in details. The present calculations 
rely strongly on the new developed 4th-order code. An explicit 
comparison between the two codes is shown in Appendix A. 
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FIG. 2. Complete I = m = 2 gravitational (Zerilli) waveform corresponding to the dynamics depicted in Fig. 1. The waveform 
is extracted at rf'^'/M = 1000. 



plunge phase, like the computation of the recoil velocity 
imparted to the center of mass of the system due to the 
linear momentum carried away by GWs [73]. As we will 
see in Sec. IV, high-order multipoles are, in fact, needed 
to obtain an accurate result. An analysis of the relative 
importance of the different multipoles based on energy 
considerations is the subject of Sec. Ill A. 

As for the instantaneous GW frequency, the right- 
panels of Fig. 3 show the same kind of behavior for 
each multipolc: Muogm is approximately equal to mil 
during the inspiral, to grow abruptly during the nona- 
diabatic plunge phase until it saturates at the ringdown 
frequency (indicated by dashed lines in the plot). As al- 
ready pointed out in Ref. [31, 32] the oscillation pattern 
that is clearly visible for some multipoles is due to the 
contemporary (but asymmetric) excitation of the posi- 
tive and negative frequency QNMs of the black hole. We 
shall give details on this phenomenon in Sec. IIIB. 

A. quasiuniversal plunge 

In this section we discuss in quantitative terms the rel- 
ative contribution of each multipole during the plunge, 
merger and ringdown phase. The analysis is based on 
the energy and angular momentum computed from the 
emitted GW. While these quantities represent a "synthe- 
sis" of the information we need, their computation and 
interpretation have some subtle points that are discussed 
below. 

For a (adiabatic) sequence of circular orbits, this in- 
formation was originally obtained in Cutler et al. [38]; 
for the radial plunge of a particle initially at rest at 
infinity, the classical work of Davis, Ruffini, Press and 
Price [60] found that about the 90% of the total energy 
is quadrupole {^ = 2) radiation, and about the 8% is 
octupolc (€ = 3) radiation. Concerning the transition 



from quasicircular inspiral to plunge, Ref. [31] performed 
a (preliminary) calculation of the total energy and angu- 
lar momentum losses during a "plunge" phase (that was 
defined by the condition r < 5.9865M, with = 0.01) 
followed by merger and ringdown, computing all the mul- 
tipolar contributions up to € = 4 (see Table 2 in [31]). 

We will follow up and improve the calculation of 
Ref. [31]. Let us first point out some conceptual diffi- 
culties. As a matter of fact, any kind of computation of 
the losses during the transition from inspiral to plunge in 
our setup will depend both on the value of v, and on the 
initial time from which one starts the integration of the 
fiuxes (for instance on the time when one defines the be- 
ginning of the "plunge" phase). It follows that, if robust 
and meaningful results arc desired, the calculation has to 
be focused on the part of the waveforms that is quasiu- 
niversal (i.e., with negligible dependence on v). As was 
pointed out in [31], the quasiuniversal behavior reached 
in the — > limit is linked to the quasigeodesic charac- 
ter of the plunge motion, which approaches the geodesic 
which starts from the LSO in the infinite past with zero 
radial velocity. 

In this respect, let us recall that, as shown in Ref. [15], 
the transition from the adiabatic inspiral to the nonadi- 
abatic plunge is not sharp, but rather blurred, namely it 
occurs in a radial domain around the LSO which scaled 
with I' as r — 6M ^ aMv^^^^, with the radial velocity 
scaling as Vr ^ —f3iy^^^. In practical terms, this means 
that the quasiuniversal, quasigeodesic plunge docs not 
really start at = 6M, but at about ro/A/ ~ 6 — av^l^ . 
In Ref. [15], using a 2.5PN Pade resummed radiation re- 
action, the coefficients a and /3 were determined to be 
ck2.5PN = 1-89 and /32.5PN = —0.072. However, since 
our setup is based on the 5PN resummed radiation re- 
action force, we do not expect those numbers to remain 
unchanged, so that they do not represent for us a re- 
liable estimate to extract the part of the waveforms we 
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FIG. 3. Multipolar structure of the waveform. The left panels exhibit the moduli; the right panels the instantaneous gravi- 
tational wave frequencies for some representative multipoles. Note the oscillation pattern during ringdown (especially in the 
£ = 2, m = 1 modulus and frequency) due to the interference between positve and negative frequency QNMs. The waveform 
refer to the v = 10^'^ mass ratio. 



are interested in. Taking a pragmatical approach, we can 
determine this quasiuniversal region by contrasting our 
simulations at different i/, so to see when the dependence 
on I' is sufficiently "small" (say at 1% level in the energy 
and angular momentum losses, see below). 

Figure 4, displays the "convergence" to the £ = m = 2 
modulus (upper panel) and frequency (lower panel) for 
the three values of v. For convenience the waveforms 
have been time-shifted so that the maxima of the wave- 
form mudulus (located at Wmax = in the figure) coin- 
cide. The plot clearly shows that the late-time part of the 
waveform has a converging trend to some "universal" pat- 
tern that progressively approximates the "exact" v ^ 
case. Note that, at the visual level, amplitudes and fre- 
quencies for 1/ = {10~^ 10~^} look barely distinguishable 
during the late part of the plunge, entailing a very weak 
dependence on the properties of radiation reaction. From 



this analysis we can assume a quasiuniversal and quasi- 
geodesic plunge starting at about u — Umax = — 50 (ver- 
tical dashed line), which corresponds to M0J22 — 0.167, 
which is about 1.23 x (2nLSo) (for reference, we indicate 
with a horizontal line the 2f2Lso frequency in the lower 
panel of the figure) ^ . We integrate the multipolar energy 
and angular momentum fiuxes from u — Umax = — 50Af 
onwards and sum over all the multipoles up to = 8. 
The outcome of this computation is listed in Table I 
for the V = {10~^, 10"'^}. Note that the agreement of 
these numbers at the level of 1% is a good indication of 



Note that the radial separation that corresponds to m — Umax = 
-50 is r ~ 5.2M [more precisely, r ~ 5.199M (5.228Af) for u = 
i.e., we have a 13% difference with the value, 
5.88 obtained using the former EOB analysis with a = 1.89. 
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(M-«max)/M 

FIG. 4. "Convergence" of the waveform when v ^ 0. Re- 
tarded times have been shifted so that the zero coincides with 
the maximum of the waveform modulus |*I'22'| for each value 
of V. The horizontal dashed line indicates the adiabatic LSO 
frequency. The vertical dashed line conventionally identifies 
the beginning of an approximately quasiuniversal and quasi- 
geodesic plunge phase. 



the quasigeodesic character of the dynamics behind the 
part of the waveform that we have selected. The numer- 
ical information of Table I is completed by Tables IX-X 
in Appendix B, were we list the values and the relative 
weight of each partial multipolar contribution. Coming 
thus to the main conclusion of this analysis, it turns out 
that the £ = m = 2 multipole contributes to the total 
energy (angular momentum) for about the 58% (62%), 
the £ = m = 3 for about the 20% (20%) , the ^ = m = 4 
for about the 8% (7.6%) and the ^ = m = 5 for the 3.5% 
(3.3%). For what concerns the odd-parity multipole, the 
dominant one, £ = 2, m = 1, contributes to 4.3% of the 
total energy and 2.3% of the total angular momentum. 
We address again the reader to Appendix B for the fully 
precise quantitative information. 



B. Ringdown 

Let us focus now on the analysis of the waveform dur- 
ing pure ringdown only. Our main aim here is to ex- 
tract quantitative information from the oscillations that 
are apparent in the gravitational wave frequency (and 
modulus) during ringdown (sec Fig. 3). As explained 
in Sec. IIIB of Rcf. [32], the physical interpretation of 
this phenomenon is clear, namely it is due to an asym- 
metric excitation of the positive and negative QNM fre- 
quencies of the black hole triggered by the "sign" of 
the particle motion (clockwise or counterclockwise). The 
modes that have the same sign of mfl are the dominant 



TABLE I. Total energy and angular momentum emitted 
during the quasiuniversal, quasigeodesic plunge phase, the 
merger and ringdown (it is defined by the condition Mw22 ^ 
0.167, see Fig. 4 ). 



V 








0.47688 


3.48918 


10-'* 


0.47097 


3.44271 



ones, while the others with opposite sign are less excited 
(smaller amplitude). Since QNMs are basically excited 
by a resonance mechanism, their strength (amplitude) 
for a given multipole (€, m) depends on their "distance" 
to the critical (real) exciting frequency ml^max of the 
source, where rJmax indicates the maximum of the or- 
bital frequency. In our setup, the particle is inspiralling 
counterclockwise (i.e., > 0), therefore the positive fre- 
quency QNMs are more excited than the negative fre- 
quency ones. The amount of (relative) excitation will de- 
pend on m. Such QNM "interference" phenomenon was 
noted and explained already in Refs. [31, 32], although 
no quantitative information was actually extracted from 
the numerical data. We perform here this quantitative 
analysis. 

The waveform during the ringdown has the structure 

^'iJ:' - E ct^n^-"'-' + E ^w--^-*, (18) 

n n 

were we use the notation of Rcfs. [32, 73], and denote 
the QNM complex frequencies with cr^ = agn ± Wfn and 
^tran '^^^ Corresponding complex amplitudes. For each 
value of I, n indicates the order of the mode, its 
inverse damping time and a;^„ its frequency. For example, 
defining a£„i„e"'*"" = CJrnnl^tmn^ ^^e presence of 
only one QNM (e.g., the fundamental one, ti = 0) the 
instantaneous frequency computed from Eq. (18) reads 

1 + "-IraO + 2a«mO COS (2wmt + d hno) ' 

This simple formula illustrates that, if the two modes 
are equally excited (a^o — 1) then there is a destruc- 
tive interference and the instantaneous frequency is zero; 
on the contrary, if one mode (say the positive one) is 
more excited than the other, the instantaneous frequency 
oscillates around a constant value that asymptotically 
tends to a;^^ when — 0. In general, one can use a 

more sofisticated version of Eq. (19), that includes vari- 
ous overtones for a given multipolar order, as a template 
to fit the instantaneous GW frequency and to measure 
the various aimn and 'dimn during the ringdown. For sim- 
plicity, we concentrate here only on the measure of agmo , 
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TABLE II. Fit of QNM interference with Eq. (19) for a 
representative sample of multipoles. The numbers refer to 
V = 10~^. The MbJtQ column lists the values of the funda- 
mental QNMs frequencies gathered from the literature [61-64] 
(see also Ref. [65] for a recent review and for highly accurate 
computations). By contrast, the primed values are obtained 
from our numerical data by fitting the ringdown frequency for 
both aimo and Mujio. Note the good consistency between the 
two methods. 



I 


m 










MoJio 


Muj'io 


2 


1 


7.2672x10" 


-2 


7.2678x10 


-2 


0.37367 


0.37369 


2 


2 


4.8476x10" 


-3 


4.848x10"- 


■3 


0.37367 


0.37361 


3 


1 


9.3403x10" 


-2 


9.3403x10" 


~2 


0.59944 


0.59944 


3 


2 


8.008x10" 


3 


8.011 X 10 


-3 


0.59944 


0.59936 


3 


3 


5.5471x10" 


-4 


5.5477x10" 


-4 


0.59944 


0.59943 


4 


1 


9.1560x10" 


-2 


9.1559x10 


-2 


0.80917 


0.80918 


4 


2 


9.1433x10" 


-3 


9.1435x10 


-3 


0.80917 


0.80917 


4 


3 


9.0473x10" 


-4 


9.0475x10" 


-4 


0.80917 


0.80917 


4 


4 


6.382x10" 


5 


6.379x10" 


-5 


0.80917 


0.80918 



and we use directly Eq. (19). To perform such a fit"* (with 
a least-square method) we consider only the part of the 
ringdown that is dominated by the fundamental (least- 
damped) QNM; i.e., the "plateau of oscillations" approx- 
imately starting at u/M = 4340 in the right-panels of 
Fig. 3. The fundamental frequency n = has been used 
as given input, and we fit for the amplitude ratio aimo 
and relative phase 9emO- The outcome of the fit for some 
multipoles is exhibited in Table IL Note that in the third 
and fourth column we list also the values that one obtains 
by fitting also for the frequency w^g. We obtain perfectly 
consistent results. Note that for the multipole £ = 8 
we were obliged to compute the frequency only in this 
way, since we could not find this number in the results 
of [65]: we obtain the value Muj^q = 1.60619. The table 
quantifies that the strongest interference pattern, that 
always occurs for m = 1 (for any £), corresponds to a 
relative contribution of the negative frequency mode of 
the order of about 9%. This trend remains true for all 
values of i. For example, we have agio = 9.54 x 10~^ and 
0710 = 9.48 x 10~^. Note finally that the presence of the 
negative mode for the £ ~ 2, in ~ 1 shows up also in the 
corresponding modulus 1^*21 with the characteristic 
oscillating pattern superposed to the exponential decay 
(see top- left panel of Fig. 3). [See also Ref. [66] for an an- 
alytical treatment of the ringdown excitation amplitudes 
during the plunge]. 



For this particular investigation we use i/ = 10-2 data. The 
reason for this choice is that, in our grid setup, the waveforms 
are practically causally disconnected by the boundaries and we 
have a longer and cleaner ringdown than in the other two cases. 



IV. GRAVITATIONAL RECOIL 

Let us now come to the computation of the gravita- 
tional recoil, or "kick", imparted to the system due to 
the anisotropic emission of gravitational radiation. The 
calculation of these kicks in general relativity has been 
carried out in a variety of ways, that before 2005 relied 
mainly on analytical and semianalytical techniques. In 
particular, let us mention that, after the pioneering cal- 
culation of Fitchett [74] and Fitchett and Detweiler [75] , 
earlier estimates included a perturbative calculation [76] , 
a close-limit calculation [77] and a post-Newtonian cal- 
culation valid during the inspiral phase only [78]. This 
latter calculation has been recently improved by bring- 
ing together post-Newtonian theory and the close limit 
approximation, yielding close agreement with purely NR 
results [79]. In addition, a first attempt to compute the 
final kick within the BOB approach [73] yielded the an- 
alytical understanding (before any numerical result were 
available) of the qualitative behavior of the kick veloc- 
ity (notably the so-called "antikick" phase), as driven 
by the intrinsically nonadiabatic character of the plunge 
phase. Such preliminary EOB calculation was then im- 
proved in [80] , which included also inputs from NR simu- 
lations. On the numerical side, after the pioneering com- 
putation of Baker et al. [81], there has been a plethora 
of computations of the kick from spinning black holes 
binaries, focusing in particular on the so-called super- 
kick configurations. By contrast, for the nonspinining 
case, Rcfs. [12, 82] represent to date the largest span 
of mass ratios for which the final kick velocity is known 
(see also Ref. [72] for the nonprecessing, equal-mass spin- 
ning case). In addition, the use of semianalytical models 
prompted a deeper understanding of the structure of the 
gravitational recoil as computed in NR simulations [83]. 
However, despite all these numerical efforts, to date there 
arc no "numerical" computations of the final recoil veloc- 
ity in the f ^ limit: the only estimates rely on fits to 
NR data of the form 

„kick ^ ^j^Vl -4iy(l + Bl^) , (20) 

with the coefficient A giving the extrapolated value in the 
1^ ^ limit [12, 82]. The aim of this section is to pro- 
vide a value of A that comes from an actual (numerical) 
computation within perturbation theory. 

It is convenient to treat the kick velocity vector im- 
parted to the system by GW emission as a complex quan- 
tity, i.e. V = Vx + ivy. By integrating Eq. (15) in time 
and by changing the sign, the (complex) velocity accu- 
mulated by the system up to a certain time t is given 

by 

v = vx+ivy = -jj J ( + iJ-P) dt'. (21) 

Since in practical situations one is always dealing with 
a finite timcscries for the linear momentum flux, it is 
not possible to begin the integration from t — —oo, but 



9 




-0.03 -0.02 -0.01 0.01 



-iio/i'^=0. 001549- i 0.0010644 



-0.04 -0.03 -0.02 -0.01 



0.01 0.02 0.03 0.04 0.05 0.06 
,2 




4200 4250 4300 4350 

u/M 



FIG. 5. Parametric plot of versus Vy (for u = 10"^) ob- 
tained from Eq. (22) with vo/v'^ = (1.549 - il.0644) x 10"^ 
The analogous plot with vo = is shown in the inset. 



FIG. 6. Time- evolution of the magnitude of the recoil velocity. 
The figure shows the monotonic "multipolar" convergence to 
the final result. The plot refers to mass ratio v = 10""^. 



rather at a finite initial time to- This then amounts in 
the need of fixing some (vectorial) integration constant 
Wo that accounts for the velocity that the system has 
acquired in evolving from t = — cx) to t = to- i-G- 

v = vo-j^j\T^+iT^)dt. (22) 

As it was emphasized in Ref. [72], the proper inclusion 
of vq is crucial to get the correct (monotonic) qualitative 
and quantitative behavior of the time evolution of the 
magnitude \v\ of the recoil velocity. Typically, not only 
the final value of \v\ may be wrong of about a 10%, but 
one can also have spurious oscillations in during the 
inspiral phase if vq is not properly determined or simply 
set to zero. See in this respect Sec. IVA of Ref. [72]. 

As in Ref. [72] , the numerical determination of vq can 
be done with the help of the "hodograph", i.e., a para- 
metric plot of the instantaneous velocity vector in the 
complex velocity plane {vx,Vy). This hodograph is dis- 
played in Fig. 5 ior 1/ = 10~^. Let us focus first on 
the inset, that exhibits the outcome of the time inte- 
gration with uo = 0. Note that the center of the in- 
spiral (corresponding to the velocity accumulated dur- 
ing the quasiadiabatic inspiral phase) is displaced with 
respect to the correct value v = (0,0), corresponding 
to Wo = at t = — c». The initial w" and Vy are de- 
termined as the translational "shifts" that one needs to 
add (in both and Vy) so that the "center" of this in- 
spiral is approximately zero. The result of this opera- 
tion led (for ly = 10"^) ^ ^ i.i549 x 10"^ and 
v^jv^ = —1.0644 X 10^'^; this is displayed in the main 
panel of Fig. 5. 

This judicious choice of the integration constant is such 

that the modulus \v\ = v^. -t- of the accumulated re- 
coil velocity grows essentially monotonically in time and 



no spurious oscillations are present during the inspiral 
phase. This is emphasized by Fig. 6. In the figure, 
we show, as a solid line, the modulus of the total ac- 
cumulated kick velocity versus observer's retarded time 
(as before, waveforms are extracted at "/M = 1000). 
This "global" computation is done including in the sum 
of Eq. (15) all the partial multipolar contribution up to 
^ = 7 (which actually means considering also the inter- 
ference terms between 1 = 1 and ^ = 8 modes). To 
guide the eye, we added a vertical dashed line locating 
the maximum of |4'22'|i that approximately corresponds 
to the dynamical time when the particle crosses the light- 
ring. Note the typical shape of \v\, with a clean local 
maximum and the so-called "antikick" behavior, that is 
qualitatively identical to the corresponding curves com- 
puted (for different mass ratios) by NR simulations (see 
for example Fig. 1 of Ref. [83] for the 2:1 mass ratio case). 

In addition to the total recoil magnitude computed up 
to ^ = 7; we display on the same plot also the "partial" 
contribution, i.e. computations of \v\ where we truncate 
the sum over i in Eq. (15) at a given value I* < 7. In 
the figure we show (depicted as various type of nonsolid 
lines) the evolution of recoil with 2 < £* < 6. Note 
that each partial-^ contribution to the linear momentum 
flux has been integrated in time (with the related choice 
of integration constants) before performing the vectorial 
sum to obtain the total v. The fact that each curve nicely 
grows monotonically without spurious oscillations during 
the late-inspiral phase is a convincing indication of the 
robustness of the procedure we used to determine Vq 
by means of hodographs^. The figure highlights at a 



Note that the procedure can actually be automatized by deter- 
mining the "baricenter" of the adiabatic inspiral in the {vx,Vy) 
plane corresponding to early evolution. 
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TABLE III. Magnitude of the final and maximum kick veloc- 
ities for the three values of v considered. The last row lists 
the values extrapolated to v = Q from v — {lO"'^, lO"''} data. 



V 






10-2 


0.043234 


0.050547 


10"^ 


0.044401 


0.052058 


IQ-* 


0.044587 


0.052298 



0.0446 0.0523 



TABLE V. Fit coefficients for the final magnitude of the kick 
velocity from NR simulations as a function of Eq. (23). See 
text for discussion. 



Data 


A 


B 


Gonzalez et al. [82] 


0.04070 


-0.9883 


Gonzalez et al. [12, 82] 


0.04396 


-1.3012 


Gonzalez et al. [12, 82]+ This work 


0.04446 


-1.3482 



TABLE IV. Recoil velocities in the test-mass limit obtained 
by (extrapolating) different finite-mass results. Our "best" 
value is shown in bold. See text for explanations. 



Reference 


1/1/2 


Gonzalez et al. [82] 


0.04 


Damour and Gopakumar [73] 


[0.010, 0.035] 


Schnittman and Buonanno [80] 


[0.018, 0.041] 


Sopuerta et al. [77] 


[0.023, 0.046] 


Le Tiec, Blanchet and Will [79] 


0.032 


This work 


0.0446 



visual level the influence of high multipoles to achieve an 
accurate results. To give some meaningful numbers, if we 
consider £* = 6, we obtain jw^^j/i^^ = 0.0437 (i.e. 1.7% 
difference), while |wf"|/i^^ = 0.0424 (4.5% difference). 

The information conveyed by Fig. 6 is completed by 
Table III, where we list the final value of the modulus 
of the recoil velocity of the center of mass (as 
well as the corresponding maximum value |f™'*'^|/i^^) ob- 
tained in our setup for the three values of v that we have 
considered. The computation of the kick for the other 
values of v is procedurally identical and thus we show 
only the final numbers. The good agreement between the 
three numbers is consistent with the interpretation that 
the recoil is almost completely determined by the nona- 
diabatic plunge phase of the system (as emphasized in 
Ref. [73]), and thus it is almost unaffected by the details 
of the inspiral phase. Because of the late-plunge con- 
sistency between waveforms that we showed above for 
1/ = {lO"'^, lO""'}, we have decided to extrapolate the 
corresponding values of the kick for v = 0. The corre- 
sponding numbers are listed (in bold) in the last row of 
Table III. 

In Table IV we compare the value of the final re- 
coil with that (extrapolated to the test-mass limit) ob- 
tained from NR simulations [82] and with semianalytical 
or seminumerical predictions, like the BOB [73, 80], the 
close-limit approximation [77] (that all give a range, with 
rather large error bars) and the recent calculation of Le 
Tiec et al. [79] based on a hybrid post-Newtonian-close- 



limit calculation 

We conclude this section by discussing in more detail 
the comparison of our result with the NR-extrapolated 
value. Since the NR-extrapolated value that we list in 
Table IV was obtained using only the data of Ref. [82] 
(without the 10:1 mass ratio simulation of [12]), we have 
decided to redo the fit with all the NR data together (that 
have been kindly given to us by the Authors) . To improve 
the sensitivity of the fit when v gets small, we first fac- 
tor out the dependence in the data (i.e., we consider 
^NR^j^2^ by continuity with the test- mass result). We 
then fit the data with the function 

fiiy) = Ay/l-Av{l + Bu) . (23) 

Table V displays the results of the fit obtained using: 
the NR data of Ref. [82] (consistent with the published 
result), first row; the joined information of Refs. [12, 82], 
second row; and the NR data of [12, 82] together with the 
test-mass result calculated in this paper. Note that the 
NR fit are perfectly consistent with the test-mass value: 
in particular, our extrapolated value 0.0446 
shows an agreement of 1.5% with the value of A obtained 
from the fit to the most complete NR information (in bold 
in Table V). 

The information of the table is completed by Fig. 7, 
that displays f{v) (as a dash-dot line) obtained from the 
complete NR data of Refs. [12, 82]. Note the visual good 
agreement between this extrapolation and the test-mass 
point when — >■ 0. For contrast, we also show on the plot 
(as a dashed line) the outcome of the fit with the simple 
Newtonian-like formula {B = 0) [74]. We also tested the 
effect of adding a quadratic correction [i.e. a term Cv^ in 
the polynomial multiplying the square root in /(i^)], but 
we found that it does not really improve the description 
of the data. 



V. CONSISTENCY CHECKS 

In this section, we finally come to the discussion of 
some internal consistency checks of our approach. These 
consist in: (i) the verification of the consistency between 
the mechanical angular momentum loss (as driven by 
our analytical, resummcd radiation reaction force) and 
the actual gravitational wave energy flux computed from 
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FIG. 7. Results of the fit of NR data of Refs. [12, 82] us- 
ing Eq. (23). Note the good agreement between the NR- 
extrapolation and our test-mass result. The bottom panel 
contains the relative difference with the data. This plot cor- 
responds to the second row of Table V, without the test-mass 
point. See text for details. 
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FIG. 8. Comparison between two angular momentum losses: 
the GW flux (solid line) computed a la Regge-Wheeler-Zerilli 
including up to ^ = 8 radiation multipoles, and the mechanical 
angular momentum loss —F^p (dash line). The two vertical 
lines correspond (from left to right) to the particle crossing 
respectively the adiabatic LSO location (r = 6A/) and the 
light-ring location (r = 3M). The plot refers to v = 10""^. 



the waves (as a follow up of a similar analysis done in 
Ref. [32]); (ii) a brief analysis of the influence on the 
(quadrupolar) waveform of the higher-order z^-dependent 
EOB corrections entering the conservative and noncon- 
servative part of the relative dynamics. 

A. Angular momentum loss 

One of the results of Ref. [32] was about the compari- 
son between the mechanical angular momentum loss pro- 
vided by the resummed radiation reaction J-^p and the 
angular momentum flux computed from the multipolar 
waveform. At that time, the main focus of Ref. [32] was 
on to use of the "exact" instantaneous gravitational wave 
angular momentum flux J [see Eq. (14)], to discriminate 
between two different expression of the 2.5PN Pade re- 
summed angular momentum flux T"^'^^^ that are degen- 
erate during the adiabatic early inspiral. In addition , in 
that setup it was also possible to: (i) check consistency 
between J and —F^p during the inspiral and early plunge; 
(ii) argue that non-quasicircular corrections in the radi- 
ation reaction are present to produce a good agreement 
between the "analytical" and the exact angular momen- 
tum fluxes also during the plunge, almost up to merger 
and (iii) show that the "exact" flux is practically insen- 
sitive to (any kind of) NQC corrections. Since we are 
now using a new radiation reaction force with respect to 
Ref. [32] , it is interesting to redo the comparison between 
the Regge-Wheeler-Zerilli "exact" flux and the "analyti- 
cal" mechanical loss computed along the relative dynam- 
ics. The result of this comparison is displayed in Fig. 8, 
that is the analogous of (part of) Fig. 2 of Ref. [32]. We 



show in the figure the mechanical angular momentum 
loss (changed of sign) —T^jv versus the mechanical time 
tjM together with the instantaneous angular momentum 

flux j jv (computed from including all contribu- 

tions up to ^ = 8) versus observer's retarded time. Note 
that we did not introduce here a possible shift between 
the mechanical time t and observer's retarded time u. 
As such a shift is certainly expected to exist, our re- 
sults should be viewed as giving a lower bound on the 
agreement between T^p and J. Note the very good visual 
agreement, not only above the LSO (vertical dashed line) 
but also htlow the LSO, and actually almost during the 
entire plunge phase. In fact, the accordance between the 
two fluxes is actually visually very good almost up to the 
merger (approximately identifled by the maximum of the 
£ = m = 2 waveform, see dash-dot line in the inset)^. 

We inspect this agreement at a more quantitative level 
in Fig. 9, where we plot the (relative) difference between 
J/v and —T^jv versus (twice) the orbital frequency. The 
inset shows the relative difference from initial frequency 
to 2rimax! where fimax is the maximum of the orbital fre- 
quency. The main panel is a close-up centered around the 



® Following the reasoning line of [32], the result displayed in the 
figure is telling us that most of the non-quasi circular correc- 
tions to the waveforms (and energy flux) are already taken into 
account automatically in our resummed flux, due to the intrin- 
sic dependence on it on pr» through the Hamiltonian, so that 
one might need only to add pragmatically corrections that are 
very small in magnitude. These issues will deserve more careful 
investigations in forthcoming studies 
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FIG. 9. Difference between mechanical angular momentum 
loss and GW energy flux shown versus twice the orbital fre- 
quency 2f2. The vertical line locates the adiabatic LSO fre- 
quency. The main panel focuses on the inspiral phase, while 
the inset shows the full range until the 2f2max, where Slmax 
indicates the maximum of orbital frequency. 



LSO frequency. Note that the relative difference is of the 
order of 10~^ during the late inspiral and the plunge, in- 
creasing at about only a 10% just before merger. We have 
performed the same analysis for v = 10^^ and v = 10~^, 
obtaining similar results. This is an indication that we 
have reached the limit of accuracy of our resummation 
procedure, limit that evidently is more apparent during 
the late part of the plunge. It is however remarkable that 
the fractional different is so small, confirming the valid- 
ity of the improved p-resummation of Ref. [41]. In this 
respect, we mention in passing that this fractional differ- 
ence can be made even smaller by further Pade resum- 
ming the residual amplitude corrections pim hi a proper 
way. This route was explored in Ref. [41] for the P22 am- 
plitude, yielding indeed better agreement with the "ex- 
act" circularized waveform amplitude. A more detailed 
analysis of these delicate issues lies out of the purpose of 
this paper, but will be investigated in future work. 



B. Influence of dynamical "self- force" //-dependent 
effects on the waveforms. 

In the work that we have presented so far we have in- 
cluded in the relative dynamics only the leading order 
part of the radiation reaction force, namely the one pro- 
portional to v. This allowed us to compute, consistently 
as shown above, Regge- Wheeler- Zerilli- type waveforms. 
In doing so we have neglected all the Hmte-v effects that 
are important in the (complete) EOB description of the 
two-body problem, that is: (i) i/-dcpendcnt corrections to 
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TABLE VI. Accumulated phase difference (computed from 
oji = 0.10799 up to CJ2 = 2r2LSO = 0.13608 ) [in radi- 
ans] between £ = m = 2 EOB waveforms. A^^^q'^^'^ is 
the phase difference accumulated between the EOBspn and 
the EOBtestmass insplunge waveforms, while A(^l°q is the 
phase difference accumulated between the EOBipn and the 
EOBtestmass insplunge waveforms. See text for more precise 
explanations. 



V 


A0L°o=™ [rad] 


A0L°o^™ [rad] 


10"^ 


3.2 


0.40 


10"=* 


3.8 


0.43 


10-^ 


4.1 


0.44 



the conservative part of the dynamics'" and (ii) higher or- 
der V dependent corrections in the nonconservative part 
of the dynamics, i.e. corrections entering in the definition 
of the angular momentum flux T^p. 

In this section we want to quantify the effects en- 
tailed by these corrections on our result. To do so, we 
switch on the "self-force" ^/-dependent corrections in the 
Hamiltonian and in the flux defining the complete EOB 
relative dynamics and we compute EOB waveforms for 
V = {10~^, 10^'^, 10~^}. Since this analysis aims at giv- 
ing us only a general quantitative idea of the effect of 
"self-force" corrections, we restrict ourself only to the 
computation of the £ = m = 2 "insplunge" waveform, 
without the matching to QNMs [26]. Note also that we 
neglect the non-quasi-circular corrections advocated in 
Eq. (5) of [26]. (See also Ref. [19]). 

For each value oi v = {10~^, 10""^, 10"''}, we com- 
pute three insplunge /i22 resummed waveforms with in- 
creasingly physical complexity. The first, EOBtestmass 
insplunge waveform, is obtained within the 0{y) approx- 
imation used so far; i.e., we set to zero all the v depen- 
dent EOB corrections in -ffsoB and in the normalized 
fiux, /din = !v)m{vif,;v = 0). The second, EOB5PN m- 
splunge waveform, is computed from the full EOB dy- 
namics, with the complete Heob and ^/-dependent (New- 
ton normalized) flux /DiN(Wi^;i') replaced in Eq. 7. The 
radial potential A(u; 05, ag, i^) is given by the Pade re- 
summed form of Eq. (2) of Ref. [26] and 05 and ag are 
EOB flexibility parameters that take into account 4PN 
and 5PN corrections in the conservative part of the dy- 
namics. They have been constrained by comparison with 
numerical results [26, 50]. Following [50], we use here the 
values 05 = —22.3 and ag = +252 as "best choice". The 
third, EOBiPN insplunge waveform, is obtained by keep- 
ing the same flux fBm{v^p]v) of the EOB5PN case, but 
only part of the EOB Hamiltonian. More precisely, we 



Those corrections come in both from the resummed EOB Hamil- 
tonian Heob with the double-square-root structure and from the 
EOB radial potential A(r,iy). 
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restrict the effective Hamiltonian Hcs at IPN level. This 
practically means using A(r; 0) = 1 — 2M/r and dropping 
the pf^ jr^ correction term that enters in J?eff at the 3PN 
level. See Eq. (1) in [50]. We compute the relative phase 
difference, accumulated between frequencies (wi, u;2), be- 
tween the EOBtcstmass waveform and the other two. We 
chose Wi = 0.10799, that corresponds to the initial (test- 
mass) GW frequency, and W2 = 2r2Lso — 0.13608. In- 
stead of comparing the waveforms versus time, wc found 
it convenient to do the following comparison versus fre- 
quency. For each waveform, we compute the following 
auxiliar quantity 

This quantity measures the effective number of GW cy- 
cles spent around GW frequency lo (and correspondingly 
weighs the signal-to- noise ratio [85]), and is a useful di- 
agnostics for comparing the relative phasing accuracy of 
various waveforms [84]. Then, the gravitational wave 
phase accumulated between frequencies (u^^io-i) 

is given by 

0(c^i,c^2) = / Qc^rflogw. (25) 

We can then define the relative depliasing accumulated 
between two waveforms as 

^C!r™ = / Ag50B„PNrfiog(^)^ (26) 

where A(g^°^"P« = qEOB„pn „ gEOBteetmass ^ The resuhs 
of this comparison are contained in Tabic VI. Note the 
influence of the correction due to the conservative part of 
the self force. Since this correction changes the location 
of the adiabatic r-LSO position [14], it entails a larger 
effect on the late-time portion of the binary dynamics 
and waveforms, resulting in a more consistent dephasing. 

VI. CONCLUSIONS 

We have presented a new calculation of the gravita- 
tional wave emission generated through the transition 
from adiabatic inspiral to plunge, merger and ringdown 
of a binary systems of nonspinning black holes in the ex- 
treme mass ratio limit. We have used a Regge-Wheeler- 
Zerilli perturbative approach completed by leading or- 
der EOB-based radiation reaction force. With respect to 
previous work, we have improved (i) on the numerical 
algorithm used to solve the Regge- Wheeler- Zerilli equa- 
tions and (ii) on the analytical definition of the improved 
EOB-resummed radiation reaction force. 

Our main achievements are listed below. 

1. We computed the complete multipolar waveform 
up to multipolar order I = %. We focused on the 
relative impact (at the level of energy and angular 



momentum losses) of the subdominant multipolcs 
during the part of the plunge that can be consid- 
ered quasiuniversal (and quasigeodesic) in good ap- 
proximation. We analyzed also the structure of the 
ringdown waveform at the quantitative level. In 
particular, we measured the relative amount of ex- 
citation of the fundamental QNMs with positive 
and negative frequency. We found that, for each 
value of I, the largest excitation of the negative 
modes always occurs for m = 1 and is of the order 
of 9% of the corresponding positive mode. 

2. The central numerical result of the paper is the 
computation of the gravitational recoil, or kick, im- 
parted to the center of mass of the system due to 
the anisotropic emission of gravitational waves. We 
have discussed the influence of high modes in the 
multipolar expansion of the recoil. We showed that 
one has to consider £ > 4 to have a ^ 10% accu- 
racy in the final kick. We found for the magnitude 
of the final and maximum recoil velocity the val- 
ues Iz;™'^!/!/^ ^ 0.0446 and ^ 0.0523. 
The value of the final recoil shows a remarkable 
agreement (< 2%) with the one extrapolated from 
a sample of NR simulations, \v^i\lv^ ~ 0.0439. 

3. The "improved rcsummation" for the radiation re- 
action used in this paper yields a better consistency 
agreement between mechanical angular momentum 
losses and gravitational wave energy flux than the 
previously employed Pade resummed procedure. In 
particular, we found an agreement between the an- 
gular momentum fluxes of the order of 0.1% during 
the plunge (well below the LSO), with a maximum 
disagreement of the order of 10% reached around 
the merger. This is a detailed piece of evidence 
that EOB waveforms computed via the resumma- 
tion procedure of [41] can yield accurate input for 
LISA-oriented science. 

While writing this paper, we became aware of a similar 
calculation of the final recoil by Sundararajan, Khanna 
and Hughes [86] . Their calculation is based on a different 
method to treat the transition from inspiral to plunge 
(see Refs. [87-89] and references therein). In the limiting 
case of a nonspinning binary, their results for the final 
and maximum kick are fully consistent with ours. 
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FIG. 10. Self-convergence factor computed from p-norms of 
{£ = 8,m = 8) inspiral plunge waveforms. Resolutions used 
are: Ar, = {0.05, 0.025, 0.0125}. The particle is initially at 
R{0) = 6.5M and u = 0.01. 
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Appendix A: Numerical framework, tests and 
comparison with the literature 

The numerical procedure adopted is similar to the one 
of Ref [31], but it has been improved on several aspects. 
In particular, the original code has been fully rewritten 
and optimized and a new finite-differencing algorithm to 
solve the Regge- Wheeler- Zerilli equation has been imple- 
mented. 

The Regge- Wheeler-Zerilli equations, Eqs. (10), are 
solved as a first-order-in-time second-ordcr-in-space sys- 
tem adopting the method of lines. Time advancing is 
done by means of a 4th order Runge-Kutta algorithm, 
while centered 4th order finite difi^erences are used to ap- 
proximate the space derivative. Standard Sommerfeld- 
maximally dissipative boundary conditions are adopted 
and implemented as described in [90] 

We solve the equations given in Sec. II for the par- 
ticle dynamics using a standard 4th order Runge-Kutta 
algorithm with adaptive stepsize. Then we insert the re- 
sulting position and momenta in the source terms s'f^l°^ 
using a Gaussian-function representation of 5{r^, — i?H.(i)) 
(see below). The distributional (5-function that appears 
in the source terms is approximated by a smooth function 



(5^(7"*). We use 



6{r, - R,{t)) ^ 6„{r,) 



1 



exp 



2o-2 

(Al) 

with (T > Ar, . In practice a ~ Ar, works well thanks to 
the effective averaging entailed by the fact that i?*(<) is 
not restricted to the r„ grid, but varies nearly continu- 
ously on the axis. In Ref. [31] it was already pointed 
out that, if a is sufficiently small and the resolution suf- 
ficiently high (so that the Gaussian function is resolved 
by a sufficiently high number of points) this technique is 
competitive with other approaches that prefer a mathe- 
matically more rigorous treatment of the J-function [67- 
69] (see in this respect Table 1 in Ref. [31]). Since in 
this paper we use a different numerical method to solve 
Eqs. (10), we have performed exnovo all the accuracy 
tests for circular orbits and radial plunge that were for- 
merly discussed in [31]. 

Self-convergence tests showed the correct convergence 
rate both in norm, see Fig. 10 for an example, and point- 
wise, both with and without the particle source. In the 
latter case however results were not satisfactory if the 
Gaussian in the source was not enough resolved. We 
found the correct convergence rate using, for the lowest 
resolution, a Gaussian width of cr > 3Ar*, optimal re- 
sults were obtained with a ~ lOAr*, while smaller values 
gives experimental rate around 2nd (a = Ar*) and 3rd 
order (cr = 2Ar*). Together with the physical require- 
ment a <^ M and the necessity of extracting waveforms 
at large radii, this fact poses some limits on the resolu- 
tion to be used and on the minimal computational time 
necessary for the simulations. As expected, no spurious 
oscillations were found in the Regge- Wheeler- Zerilli solu- 
tion in the region across the (smoothed) delta function. 

A direct comparison with the old code. Fig. 11, shows 
clearly that the numerical improvements lead to quanti- 
tative better results. The amplitude of '^22 is computed 
for two different resolution Ar* = 0.5 and 0.25 using the 
old code based on the Lax-Wendroff scheme and with 
Ar* = 0.5 with the new scheme. 

The differences are due only to the numerical scheme 
employed for the wave equation, since also the new radia- 
tion reaction has been used in the old code. As expected, 
the new numerical scheme shows a faster convergence and 
strongly suppresses the spurious oscillations coming from 
the boundaries; note in this respect the small "bumps" 
at t ^ 400 that arc present only in the data computed 
with the old code. 

To validate the physical results of the code at a more 
quantitative level we performed "standard" comparisons 
with the literature considering circular orbits and radial 
plunge, following the line of [31]. 

The energy and angular momentum fluxes computed 
from the waveforms generated by a particle on a circular 
orbit of radius vq ~ 7.9456M are displayed in Table Vll 
and Table Vlll. The numbers are compared with those 
present in the literature, showing very good agreement. 
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TABLE VII. Gravitational wave (multipolar) energy flux Etml^i^ of a particle on a circular orbit of radius ro = 7.9456M 
Comparison between our results and those present in the literature. Our waveforms are extracted at r, = lOOOAf. 



i m 


This work 


Ref. [69] 




Diff.[%] 


Ref [91] 




Diff.[%] 


Ref. [92] 




Diff.[%] 


2 1 


8.1733x10"^ 


8.1662x10" 


-7 


0.086 


8.1623x10" 


-1 


0.134 


8.1633x10" 


-7 


0.122 


2 2 


1.7069x10"* 


1.7064x10" 


-4 


0.029 


1.7051x10" 


-4 


0.105 


1.7063x10" 


-4 


0.035 


3 1 


2.1785x10"^* 


2.1732x10" 


-9 


0.242 


2.1741x10" 


-9 


0.200 


2.1731x10" 


-9 


0.246 


3 2 


2.5218x10"^ 


2.5204x10' 


-7 


0.057 


2.5164x10" 


-7 


0.217 


2.5199x10" 


-7 


0.077 


3 3 


2.5483x10"^ 


2.5475x10" 


-5 


0.031 


2.5432x10" 


-5 


0.201 


2.5471x10" 


-5 


0.047 


4 1 


8.3699x10"^^ 


8.4055x10" 


13 


0.423 


8.3507x10" 


13 


0.230 


8.3956x10" 


13 


0.306 


42 


2.5125x10"^ 


2.5099x10" 


-9 


0.103 


2.4986x10" 


-9 


0.556 


2.5091x10" 


-9 


0.135 


4 3 


5.7792x10"** 


5.7765x10" 


-8 


0.046 


5.7464x10" 


-8 


0.571 


5.7751x10" 


-8 


0.071 


44 


4.7283x10"^ 


4.7270x10" 


-6 


0.027 


4.7080x10" 


-6 


0.430 


4.7256x10" 


-6 


0.056 


5 1 


1.2904x10-1'' 


1.2607x10" 


15 


2.357 


1.2544x10" 


15 


2.871 


1.2594x10" 


15 


2.462 


5 2 


2.7874x10"^^ 


2.7909x10" 


12 


0.126 


2.7587x10" 


12 


1.040 


2.7896x10" 


12 


0.080 


5 3 


1.0946x10"^ 


1.0936x10" 


-9 


0.095 


1.0830x10" 


-9 


1.074 


1.0933x10" 


-9 


0.122 


5 4 


1.2334x10"** 


1.2329x10" 


-8 


0.038 


1.2193x10" 


-8 


1.154 


1.2324x10" 


-8 


0.078 


5 5 


9.4630xl0~^ 


9.4616x10" 


-7 


0.014 


9.3835x10" 


-7 


0.847 


9.4563x10" 


-7 


0.071 




350 400 450 500 550 600 650 700 
t/M 

FIG. 11. Comparison between the old and the new code. The 
plot shows the amplitude of the (2, 2) mode of the Regge- 
Wheeler-Zerilli as computed with the old code (at resolutions 
Ar* — 0.5 and 0.25) based on the Lax-Wendroff scheme and 
with the new code (at resolution Ar, = 0.5). Other parame- 
ters of the runs are jj, — 0.001, -R(O) = 6.5 and Cca = 0.9. 




-0.3 - 
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u/{2M) 

FIG. 12. Waveform emitted by a particle plunging radially 
into the black hole (along the z-axis) from r = lOM. The 
waveform is extracted at r* — lOOOM. 



The fractional differences are always well below the 1% 
except for the multipole (5, 1) (2% in the energy flux) 
whose absolute value is the smallest. Notice that, dif- 
ferently from Ref. [31], the accuracy is maintained also 
for high multipoles. Wc also computed multipolcs for 
£ > 6, although we did not report them here since cor- 
responding data to compare with arc not present in the 
literature. 

For what concerns the radial infall. Fig. 12 displays the 
£ = 2, m = waveform generated by a particle plunging 
into the black hole radially along the z-axis. The parti- 



cle has zero initial velocity and starts at r = lOM. We 
specify conformally flat initial data according to the pro- 
cedure described at the end of Sec. 4 of [31]. Note that 
^'20^ has been multiplied by a factor 2 to facilitate the 
(very satisfactory) comparison with the top-right panel 
of Fig. 4 in [68] and top- left panel of Fig. 6 in [67] . 

An additional test involved the dependence of the re- 
sults on the Gaussian ampitude a. Focusing on fj, = 0.01 
and inspiral plunge simulations with i?(0) = 6.5 M, we 
experimented with the values a = {0.05,0.1,0.5,1} M 
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TABLE VIII. Gravitational wave (multipolar) angolar momentum flux Jim/y? of a particle on a circular orbit of radius 
To = 7.9456M. Comparison between our results and those present in the literature. Our waveforms are extracted at r» = lOOOM. 



I m 


This work 


Ref. [69] 




Difr.[%] 


Ref. [91] 




Diff.[%] 


Ref. [92] 




Difr.[%] 


2 1 


1.8305x10"^ 


1.8289x10" 


-5 


0.090 


1.8270x10' 


-5 


0.194 


1.8283x10' 


-5 


0.122 


2 2 


3.8229x10"^ 


3.8219x10' 


-3 


0.027 


3.8164x10' 


-3 


0.171 


3.8215x10' 


-3 


0.037 


3 1 


4.8790x10"** 


4.8675x10' 


-8 


0.237 


4.8684x10' 


-8 


0.219 


4.8670x10' 


-8 


0.247 


3 2 


5.6481x10"'' 


5.6450x10" 


-6 


0.055 


5.6262x10' 


-6 


0.389 


5.6439x10' 


-6 


0.075 


3 3 


5.7074x10"* 


5.7057x10' 


-4 


0.030 


5.6878x10' 


-4 


0.345 


5.7048x10' 


-4 


0.046 


4 1 


1.8518x10"" 


1.8825x10" 


11 


1.633 


1.8692x10" 


-11 


0.933 


1.8803x10" 


-11 


1.517 


4 2 


5.6272x10"* 


5.6215x10' 


-8 


0.102 


5.5926x10' 


-8 


0.619 


5.6195x10' 


-8 


0.138 


4 3 


1.2944x10"*^ 


1.2937x10' 


-6 


0.052 


1.2933x10' 


-6 


0.083 


1.2934x10' 


-6 


0.075 


4 4 


1.0590x10"* 


1.0586x10' 


-4 


0.037 


1.0518x10' 


-4 


0.684 


1.0584x10' 


-4 


0.056 


5 1 


2.8558x10"" 


2.8237x10" 


-14 


1.138 


2.8090x10" 


-14 


1.667 


2.8206x10" 


-14 


1.249 


5 2 


6.2386x10"" 


6.2509x10" 


11 


0.197 


6.1679x10" 


-11 


1.146 


6.2479x10" 


-11 


0.149 


5 3 


2.4517x10"* 


2.4494x10' 


-8 


0.093 


2.4227x10' 


-8 


1.196 


2.4486x10' 


-8 


0.125 


5 4 


2.7625x10"'' 


2.7613x10' 


-7 


0.042 


2.7114x10' 


-7 


1.883 


2.7603x10' 


-7 


0.078 


5 5 


2.1194x10"^ 


2.1190x10' 


~5 


0.020 


2.0933x10' 


- 5 


1.248 


2.1179x10' 


-5 


0.072 



using typical resolutions. The extreme value a = M 
(see discussion above) gave reliable waveforms, while 
spurious modulations due to the extended source were 
clearly evident for a > M . Relative differences between 
(2, 2)-waveforms, taking a = 0.05 M as reference value, 
were of the order of {10~^, 10~^, 10"''} in the ampli- 
tudes respectively for tr = {0.01,0.5,1} M, and 1 order 
of magnitude less for the phase. For the (8,8) multipole 
(worst case) they were {10~^, 10"^, 10"'^} always for cr = 
{0.01,0.5,1} M. The relative differences on the energy 
flux multipoles were of the order of {0.0002, 0.007, 0.03} 
in the (2, 2) case and {0.005, 0.2, 0.9} in the (8, 8) case. 
While a > 0.5 M docs not give satisfactory results, differ- 
ences between a = 0.05 and ct = 0.1 are quite small, giv- 
ing the same results up to 0.5% for the energy contribute. 
In the simulations of the paper we used a = 0.04M. 

We also tested the use of the alternative source term 
(mathematically equivalent for a distributional source) 
given by Eq (23) of [31]. Consistently with the analysis 
of Ref. [31] (see their Fig. 5), for insplunge waveforms 
we found a relative difference smaller than 10"^ both in 
amplitude and phase in the £ = m = 2 mode. 

We finally mention that in the case of long simula- 
tions (i.e., v = 10"^, ro = 7M) our boundary conditions 
are not fully satisfactory. In fact, in this case the wave- 
forms might be (slightly) contaminated by small reflec- 
tions from boundaries, especially for high multipoles. A 
solution to this problem is discussed in [93] (and refer- 



ences therein), in which, basically, the free-data in the 
Sommcrfeld condition (in our case they are set to zero) 
are specified as an integral convolution between a time- 
domain boundary kernel and the solution. The method 
provides an exact radiative outer boundary condition for 
the wave equations. We mention that an alternative 
approach is represented by solving the Regge-Wheeler- 
Zerilli equation on matched hyperboloidal foliations as 
described in [94] (and references therein) . The appealing 
feature is the fact that it has the double advantage to 
not need boundary conditions (no incoming modes) and 
to allow extraction exactly at null infinity. 



Appendix B: Partial losses during the plunge phase 

Let us list here some useful numerical information re- 
lated to the discussion of Sec. Ill A; i.c, the energy and 
angular momentum emitted during the quasiuniversal 
part of the late-plunge, merger and ringdown. The pure 
numbers, up to multuipolar order ^ = 8, are given in 
Table IX. The relative percentage (with respect to the 
"total" energy -Etot and angular momentum Jtot) are 
given in the following Table X. Note here that by total we 
indicate the sum over multipoles (€, m), with 2 < £ < 8 
and \m\ < I. 



[1] F. Pretorius, Phys. Rev. Lett. 95, 121101 (2005) [2] M. Campanelli, C. O. Lousto, P. Marronetti and Y. Zlo- 
[arXiv:gr-qc/0507014]. chower, Phys. Rev. Lett. 96, 111101 (2006) [arXiv:gr- 

qc/0511048]. 
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TABLE IX. Multipolar contributions to total energy and angular momentum emitted during the quasiuniversal part of the 
plunge phase, the merger and the ringdown. 





Multipolc 


/i = 0.001 




/i = 0.0001 




I 


m 




Jiml 




Jlml 


2 





9.739x10-"* 





9.739x10"'' 







1 


2.023 xlO"^ 


0.787 xlO"^ 


2.035 xlO"^ 


0.793 xlO"^ 




2 


2.733 xlO^* 


2.150 


2.769 xlO"* 


2.179 


3 





3.320 xlO"^ 





3.330x10"^ 







1 


0.540 xlO"^ 


1.196 xlO"^ 


0.543 xlO^^* 


1.201 xlQ-^* 




2 


0.782 xlO"^ 


0.370 xlO-i 


0.788 xlO"^ 


0.373 xlO"^ 




3 


0.929 xlO"^ 


0.678 


0.941 xlO"^ 


0.687 


4 





1.518 xlO"** 





1.523x10"'* 







1 


2.226 xlO"^ 


0.333 xlO"^ 


2.234 xlO"^ 


0.334 xlO"* 




2 


3.023 xlO""* 


0.998 xlO"^ 


3.039 xlO"'' 


1.003 xlO"^ 




3 


0.330 xlO-2 


1.657 xlO^^ 


0.332 xlO"^ 


1.674 xlO"^ 




4 


0.377 xlO"* 


2.627 xlO~^ 


0.382 xlO"* 


2.663 xlO"* 


5 





7.932 xlO-* 





7.957 xlO"* 







1 


1.042 xlO"® 


1.205 xlO^^ 


1.046 xlO"" 


1.209 xlO"^ 




2 


1.341 xlO"^ 


0.323 xlO"'' 


1.347 xlO"'^ 


0.325 xlO-* 




3 


1.659 xlO"* 


0.642 xlO"^ 


1.670 xlO"* 


0.647 xlO"^ 




4 


1.478 xlO"^ 


0.762 xlO"^ 


1.493 xlO"^ 


0.770 xlO-2 




5 


1.683 xlO"^ 


1.129 xlO"^ 


1.705 xlO"^ 


1.145 xlO"^ 


6 





4.182x10"^ 





4.197x10"^* 







1 


0.551 xlO"^ 


0.516 xlO"'' 


0.552 xlO"'^ 


0.518 xlO"'^ 




2 


0.673 xlO"*^ 


1.312 xlO"*^ 


0.676 xlO-** 


1.317 xlO"** 




3 


0.782 xlO"^ 


2.340 xlO-'^ 


0.786 xlO"^ 


2.352 xlO"'^ 




4 


0.904 xlO^"* 


0.379 xlO"^ 


0.910 xlO"'' 


0.382 xlO^^* 




5 


0.694 xlO"^ 


0.361 xlO"^ 


0.701 xlO"^ 


0.365 xlO~^ 




6 


0.799 xlO"^ 


0.518 xlO"^ 


0.809 xlO"^ 


0.526 xlO"* 


7 





2.416x10"^° 





2.424 







1 


2.948 xlO"^ 


2.337 xlO"^ 


2.959 xlO"^ 


2.345 xlO"^ 




2 


0.360 xlO"^ 


0.586 xlO^'' 


0.362 xlO"'^ 


0.588 xlO"'' 




3 


0.423 xlO"*^ 


1.058 xlO"'' 


0.425 xlO"" 


1.063 xlO"*^ 




4 


0.448 xlO"^ 


1.516 xlO^^ 


0.451 xlO"^ 


1.526 xlO"^ 




5 


0.492 xlO"* 


2.168 xlO^"* 


0.496 xlO"" 


2.186 xlO"" 




6 


0.337 xlO"'* 


1.760 xlO"^ 


0.341 xlO"^ 


1.781 xlO"^ 




7 


0.396 xlO-2 


2.497 xlO~^ 


0.401 xlO"^ 


2.533 xlO"^ 


8 





1.359 xlO"" 





1.364 xlO"" 







1 


1.699 xlO"^° 


1.163 xlO-^" 


1.704x10"^° 


1.167x10"'° 




2 


1.986 xlO"^ 


2.792 xlO"" 


1.994 xlO"^ 


2.802 xlO"^ 




3 


2.303 xlO"* 


0.491 xlO"'' 


2.313 xlO"* 


0.493 xlO"'' 




4 


2.604 xlO"'^ 


0.756 xlO"'' 


2.618 xlO"^ 


0.760 xlO"*^ 




5 


2.539 xlO"'' 


0.931 xlO"^ 


2.558 xlO"" 


0.938 xlO"'^ 
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TABLE X. Relative contribution of each multipole to the total energy and angular momentum emitted during the quasiuniversal 
part of the plunge phase, the merger and the ringdown. 





Multipole 




0.001 




= 0.0001 
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-E£m/-ETOT[%] 
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0.2050 
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2.2868 
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2.2730 
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62.45 
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62.46 
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7.0504 xlO"^ 





6.9844 xlO"^ 
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0.3475 xlO"^ 


1.1378 xlO"^ 


0.3443 xlO"^ 
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0.6972 


0.4799 




4 


8.006 


7.631 


8.007 


7.633 


5 





1.6843 xlO"'^ 





1.6686x10"^^ 







1 


2.2133 xlO"* 


0.3501 xlO"-* 


2.1936 xlO"" 


0.3466 xlO"" 




2 


2.8483 xlO"^ 


0.9389 xlO"^ 


2.8252 xlO"^ 


0.9305 xlO"^ 




3 


0.3523 xlO"^ 


0.1865 xlO"^ 


0.3501 xlO"i 


0.1853 xlO"i 




4 


0.3139 


0.2213 


0.3130 


0.2208 




5 


3.574 


3.280 


3.575 


3.281 


6 





8.8792x10"'^ 





8.8013x10"^^ 







1 


1.1691 xlO"^ 


1.4994 xlO"^ 


1.1585 xlO"'^ 


1.148428 xlO"® 




2 


1.4293 xlO"* 


0.3810 xlO"* 


1.4172 xlO"'' 


0.3774 xlO"" 




3 


1.6611 xlO"^ 


0.6796 xlO"^ 


1.6492 xlO"^ 


0.6742 xlO"^ 




4 


1.9185 xlO"^ 


1.1019 xlO"^ 


1.9085 xlO"^ 


1.0957 xlO"^ 




5 


1.4729 xlO"^ 


1.0489 xlO"^ 


1.4703 xlO"^ 


1.0470 xlO"^ 




6 


1.6963 


1.5060 


1.6970 


1.5069 


7 





5.1297x10"* 





5.0836x10"® 







1 


0.6260 xlO"® 


0.6787 xlO"'' 


0.6205 xlO"** 


0.6720 xlO"'' 




2 


0.7651 xlO"^ 


1.7025 xlO"*' 


0.7585 xlO""^ 


1.6861 xlO"® 




3 


0.8986 xlO"'' 


3.0728 xlO"^ 


0.8916 xlO"'' 


3.0458 xlO"'^ 




4 


0.9508 xlO"^ 


0.4404 xlO"^ 


0.9450 xlO"^ 


0.4373 xlO"^ 




5 


1.0456 xlO"^ 


0.6296 xlO"^ 


1.0411 xlO"^ 


0.6266 xlO"^ 




6 


0.7152 xlO"^ 


0.5111 xlO"^ 


0.7145 xlO"^ 


0.5106 xlO"^ 




7 


0.8405 


0.7254 


0.8412 


0.7260 


8 





2.8851x10"" 





2.8611x10"'^ 







1 


0.3608 xlO"^ 


0.3378 xlO"** 


0.3574 xlO"'' 


0.3344 xlO"** 




2 


0.4217 xlO"^ 


0.8109 xlO"'' 


0.4180 xlO"'' 


0.8031 xlO"'' 




3 


0.4889 xlO"^ 


1.4261 xlO"^ 


0.4850 xlO"'^ 


1.4132 xlO"® 




4 


0.5529 xlO"* 


2.1963 xlO"'^ 


0.5490 xlO"" 


2.1788 xlO"'' 




5 


0.5391 xlO"'* 


2.7038 xlO"-* 


0.5364 xlO"^ 


2.6880 xlO"" 




6 


0.5725 xlO"^ 


0.3553 xlO"^ 


0.5706 xlO"^ 


0.3539 xlO"^ 




7 


0.3566 xlO"^ 


0.2551 xlO"^ 


0.3566 xlO"^ 


0.2550 xlO"^ 




8 


0.4302 


0.3624 


0.4307 


0.3628 



20 



Phys. 8, 2155 (1967). 

[60] M. Davis, R. Ruffini, W. H. Press and R. H. Price, Phys. 
Rev. Lett. 27, 1466 (1971). 

[61] S. Chandrasekhar, "The mathematical theory of black 
holes," Oxford, UK: Clarendon (1992) 646 p. 

[62] S. Chandrasekhar and S. Detweiler, Proc. Roy. Soc. 
Lond. A 344, 441 (1975). 

[63] E. W. Leaver, Proc. Roy. Soc. Lond. A 402, 285 (1985). 

[64] E. Berti, V. Cardoso and C. M. Will, Phys. Rev. D 73, 
064030 (2006) [arXiv:gr-qc/0512160]. 

[65] E. Berti, V. Cardoso and A. O. Starinets, Class. 
Quant. Grav. 26, 163001 (2009) [arXiv:0905.2975 
[gr-qc]]. Note that Kerr black- hole frequen- 
cies (up to £ = 7) can be downloaded here 
http: / / wugrav.wustl.edu / people/BERTI/qnms.html 

[66] S. Hadar and B. Kol, arXiv:0911.3899 [gr-qc]. 

[67] C. O. Lousto and R. H. Price, Phys. Rev. D 56, 6439 
(1997) [arXiv:gr-qc/9705071]. 

[68] K. Martel and E. Poisson, Phys. Rev. D 66, 084001 
(2002) [arXiv:gr-qc/0107104]. 

[69] C. F. Sopuerta and P. Laguna, Phys. Rev. D 73, 044028 
(2006) [arXiv:gr-qc/0512028]. 

[70] K. S. Thorne, Rev. Mod. Phys. 52, 299 (1980). 

[71] M. Ruiz, R. Takahashi, M. Alcubierre and D. Nunez, 
Gen. Rel. Grav. 40, 2467 (2008) [arXiv:0707.4654 [gr- 
qc]]. 

[72] D. PoUney et al, Phys. Rev. D 76, 124002 (2007) 

[arXiv:0707.2559 [gr-qc]]. 
[73] T. Damour and A. Gopakumar, Phys. Rev. D 73, 124006 

(2006) [arXiv:gr-qc/0602117]. 
[74] M. J. Fitchett, Mon. Not. Roy. Astron. Soc. 203, 1042 

(1983). 

[75] M. J. Fitchett and S. Detweiler, Mon. Not. Roy. Astron. 

Soc. 211, 933 (1984). 
[76] M. Favata, S. A. Hughes and D. E. Holz, Astrophys. J. 

607, L5 (2004) [arXiv:astro-ph/0402056]. 
[77] C. F. Sopuerta, N. Yunes and P. Laguna, Phys. Rev. 

D 74, 124010 (2006) [Erratum-ibid. D 75, 069903 (2007 



ERRAT,D78,049901.2008)] [arXiv:astro-ph/0608600]. 

[78] L. Blanchet, M. S. S. Qusailah and C. M. Will, Astro- 
phys. J. 635, 508 (2005) [arXiv:astro-ph/0507692]. 

[79] A. Le Tiec, L. Blanchet and C. M. Will, Glass. Quant. 
Grav. 27, 012001 (2010) [arXiv:0910.4594 [gr-qc]]. 

[80] J. D. Schnittman and A. Buonanno, arXiv:astro- 
ph/0702641. 

[81] J. G. Baker, J. CentreUa, D. I. Choi, M. Koppitz, 

J. R. van Meter and M. C. Miller, Astrophys. J. 653, 

L93 (2006) [arXiv:astro-ph/0603204]. 
[82] J. A. Gonzalez, U. Sperhake, B. Bruegmann, M. Han- 

nam and S. Husa, Phys. Rev. Lett. 98, 091101 (2007) 

[arXiv:gr-qc/0610154]. 
[83] J. D. Schnittman et al, Phys. Rev. D 77, 044031 (2008) 

[arXiv:0707.0301 [gr-qc]]. 
[84] T. Damour, A. Nagar and M. Trias, in preparation 

(2010). 

[85] T. Damour, B. R. Iyer and B. S. Sathyaprakash, Phys. 

Rev. D 62, 084036 (2000) [arXiv:gr-qc/0001023]. 
[86] P. A. Sundararajan, G. Khanna and S. A. Hughes, 

arXiv: 1003.0485 [gr-qc]. 
[87] P. A. Sundararajan, G. Khanna and S. A. Hughes, Phys. 

Rev. D 76, 104005 (2007) [arXiv:gr-qc/0703028]. 
[88] P. A. Sundararajan, Phys. Rev. D 77, 124050 (2008) 

[arXiv:0803.4482 [gr-qc]]. 
[89] P. A. Sundararajan, G. Khanna, S. A. Hughes 

and S. Drasco, Phys. Rev. D 78, 024022 (2008) 

[arXiv:0803.0317 [gr-qc]]. 
[90] G. Galabrese and C. Gundlach, Class. Quant. Grav. 23, 

S343 (2006) [arXiv:gr-qc/0509119]. 
[91] K. Martel, Phys. Rev. D 69, 044025 (2004) [arXiv:gr- 

qc/0311017]. 

[92] E. Poisson, Phys. Rev. D 52, 5719 (1995); Phys. Rev. D 

55,7980 (1997). 
[93] S. R. Lau, J. Math. Phys. 46, 102503 (2005) [arXiv:gr- 

qc/0507140]. 

[94] A. Zenginoglu, Class. Quant. Grav. 27, 045015 (2010) 
[arXiv:0911.2450 [gr-qc]]. 



